Intrinsic plasmons in 2D Dirac materials 
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We consider theoretically, using the random phase approximation (RPA), low-energy intrinsic 
plasmons for two-dimensional (2D) systems obeying Dirac-like linear chiral dispersion with the 
chemical potential set precisely at the charge neutral Dirac point. The "intrinsic Dirac plasmon" 
energy has the characteristic dispersion in the 2D wave-vector q, but vanishes as y/T in temper- 
ature for both monolayer and bilayer graphene. The intrinsic plasmon becomes overdamped for a 
fixed q as T — >■ since the level broadening (i.e. the decay of the plasmon into electron-hole pairs 
due to Landau damping) increases as 1/VT as temperature decreases, however, the plasmon mode 
remains well-defined at any fixed T (no matter how small) as g — >■ 0. We find the intrinsic plasmon 
to be well-defined as long as g < . We give analytical results for low and high temperatures, and 
numerical RPA results for arbitrary temperatures, and consider both single-layer and double-layer 
intrinsic Dirac plasmons. We provide extensive comparison and contrast between intrinsic and ex- 
trinsic graphene plasmons, and critically discuss the prospects for experimentally observing intrinsic 
Dirac point graphene plasmons. 



I. INTRODUCTION 



Collective plasma oscillations of free carriers in doped 
or gated graphene^ (we would refer to this situation as 
"extrinsic" graphene where the chemical potential or the 
Fermi level is doped away from the Dirac point) have 
attracted considerable interest both from fundamental 
and technological perspectives^Ti^. The fundamental 
interest arises from the fact that graphene plasmons 
are apriori quantum-mechanical entities with no classi- 
cal analogs whatsoever— (since in classical physics en- 
ergy is always proportional to the square of the momen- 
tum and never has a linear dispersion as in graphene). 
This is manifested in the fact that the long wave-length 
plasma dispersion relation in extrinsic graphene goes as 
LOp = \/2/rJhvFqEp , where Ep = hvphp = Hvpy/T^n is 
the Fermi energy (i.e. the chemical potential at T = 0) 
associated with a doping carrier density of n and vp 
is the (constant) graphene velocity defining the linear 
energy dispersion, and Vs = /(nhvp) is the so-called 
graphene fine-structure constant (defining the dimension- 
less strength of Coulomb interaction with k, being the 
background lattice dielectric constant), with an "ft" ap- 
pearing explicitly in the definition of the plasma fre- 
quency [ujp = yj2e'^qvph^Jim.fK in terms of the exper- 
imentally controlled variables g, n, and k). This is in 
sharp contrast to the corresponding parabolic disper- 
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with an effective mass m where Wr, 
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) is the same classically or quantum-mechanically 

Km 

in the long wavelength [q — >■ 0) limit in any dimensional- 
ity. The technological interest arises from the consider- 
able recent progress in graphene nanoplasmonics^^— for 
prospective optoelectronic applications^"—. 

In contrast to the extensively studied extrinsic 
graphene plasmons, there has been little interest in the 
collective modes of intrinsic or undoped graphene, where 



the chemical potential sits right at the Dirac point with 
a completely filled valence band a completely empty con- 
duction band at T = 0. (Our interest here is in low- 
energy ~ meV two-dimensional collective modes and not 
very high-energy ^ 10 eV band- or so-called tt- plas- 
mons where the whole valence band charge response is 
involved^"—.) By definition, the doping carrier density 
vanishes at the Dirac point with n = (oc fc^? cx Ep)^ and 
the extrinsic graphene dispersion relation, tOp oc rt^/^, im- 
plies that no intrinsic graphene plasmons (or more gen- 
erally, Dirac plasmons) are possible. 

The above is certainly true strictly at T = where 
there can be no free carriers for Ep = 0. But, for non- 
zero temperatures, T ^ 0, the gapless nature of graphene 
leads to a thermal population of electrons (in the conduc- 
tion band) and holes (in the valence band) with equal 
density (rig = rih = n). This thermal electron- hole exci- 
tation process is known to be a power law (in fact n oc T^) 
due to the gaplessness of graphene^. Putting n oc in 
the formula for the graphene extrinsic plasma frequency, 
we conclude that there should be a finite-temperature in- 
trinsic graphene plasmon with a long-wavelength plasma 
frequency going as cOp ~ VqT. But, finite temperature 
implies that the collective mode will decay into electron- 
hole pairs even at long wavelength, and therefore, such an 
intrinsic Dirac plasmon may be ill-defined even for q 
since its decay (i.e. damping) rate (or level broadening) 
7 could exceed the mode frequency ujp making the mode 
an overdamped excitation of little interest. 

In the current work, we theoretically study intrinsic 
Dirac plasmons for both monolayer and bilayer graphene 
and for single- and double-layer systems. We obtain 
both asymptotic theoretical analytical results at long 
wavelengths and low/high temperatures, and quantita- 
tive numerical results for arbitrary wavevectors and tem- 
peratures. We use the RPA approach which should be 
well-valid for graphene plasmons at arbitrary wavevec- 
tors by virtue of its relatively small value of (< 1 
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typically). We compare and contrast the temperature de- 
pendence of intrinsic and extrinsic graphene plasmon fre- 
quency (and their Landau damping) in order to comment 
on the feasibility of the experimental observation of our 
theoretical predictions. The possible existence of high- 
temperature intrinsic graphene plasmons (with Ep ~ 0) 
is a qualitative difference between graphene and gapped 
2D semiconductor-based electron/hole systems. 



II. THEORY AND RESULTS 

Within the RPA, the collective plasmon modes of an 
electron system are given by the zeros of the complex 
dielectric function e{q,uj): 



:{q,Lj) = l-V{q)n{q,Lu)^0 



(1) 



where V{q) is the relevant bare electron-electron (i.e. 
Coulomb) interaction and Il{q,uj) is the noninteracting 
polarizability of the system. (As an aside we note that 
e{q,u;) = 1 — V{q)Il{q,uj)) is the exact expression for 
the microscopic dielectric function if Tl{q, io) is the exact 
interacting irreducible polarizability function which, of 
course, is unknown — RPA consists of replacing the ex- 
act n(g, w) by the corresponding non-interacting or bare 
polarizability function.) Equation ([Ij, as it stands, ap- 
plies for a single-component (i.e. single-layer in our case) 
system — for the double-layer case e{q,u!) should be in- 
terpreted as a matrix in the layer index with Eq. ([T|) being 
interpreted as a determinantal equation |1 — VU\ = 0. 

In general, Eq. ([!]) will have complex solutions in fre- 
quency, uj = t^p{q) + ijiq), with uip and 7 being respec- 
tively the collective mode (i.e. plasma) frequency (which 
we will often refer to as the plasmon) and its damping. If 
LUp ^ 7, the plasmon collective mode is well-defined, and 
by contrast for 7 ^ Wp, the plasmon is heavily damped 
(or even overdamped) and is not particularly relevant 
experimentally as a self-sustaining normal mode of the 
system. 

Before proceeding with the theoretical details for in- 
trinsic graphene collective modes, we write down the 
full formal expression for the graphene non-interacting 
polarizabilityii^ to be used in Eq. ([I}: 
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where A is the area of the 2D layers and the factor of 
4 arises from the valley/spin degeneracy (two each) of 
graphene. In Eq. ([2]), k' = k + q and s,s' = ±1 with 
i^s.is' (k, k') = (14- cos(?)/2 arising from the matrix el- 
ement effect associated with the chiral nature of Dirac 
fermions. The functions es_k, Cs'.k' are single-particle en- 
ergies for wavevector k, k' respectively, and /s,k, /s'.k' 
are the corresponding non-interacting Fermi distribution 
functions. We mention that explicit forms for the po- 
larizability function in monolayer and bilayer graphene 



were derived in Refs. [l| and [g^ respectively at the 
zero temperature. The finite-temperature polarizability, 
which cannot be obtained in a closed analytic form for 
arbitrary T, can be directly obtained from Eq. ^ using 
finite-temperature Fermi distribution functions or (nu- 
merically more conveniently) by using the following inte- 
gral identity to obtain the finite-T polarizability from its 
known analytic forni^ at T = 



n(g,w,^;r) = / d/i' 



n(g,c^,A^^;r = o) 

4fcBrcosh2[(/^ - ii')l2kBT\ 



(3) 



We also note that the 2D Coulomb interaction is given 
by: 



y{<i) 
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where z = for single-layer systems and \z\ = d, where 
d is the interlayer separation for double layer systems. 



A. Monolayer graphene 



0.35 

0.30 

^ 025 
H 

.J 020 

010 
005 



(a) 






K = 1 




K — h -— - " 






/ / 





(b) 




K = 1 ^ 








/ 




/ 




Y 





0.01 

AVf q/«rjT) 



0.5 1.0 1.5 

ftVf q/(A-sT) 



FIG. 1: (Color online). Plasmon dispersion of intrinsic MLG. 
(a) presents the plasmon dispersion in the high temperature 
regime (fcsT 3> hvpq)- The dashed and solid Unes correspond 
to the analytical results (given in Eq. (O) and the numerical 
results, respectively, (b) Numerical results for the plasmon 
dispersion as a function of hvFq/{kBT). The dashed and solid 
lines correspond to k = 1 and k = 5, respectively. 



We first consider the intrinsic plasmon modes in mono- 
layer graphene (MLG). We use units such that h = 
1 throughout so that frequency /energy and wavevec- 
tor/momcntum have the same units in our notations. 

For intrinsic graphene, there are no free carriers at zero 
temperature and the chemical potential fi is precisely at 
the Dirac point: ^ = 0. We note that this is true even 
for T ^ 0, i.e. fi{T) = /i(0) = Ep ^ for intrinsic 
graphene, by definition. As already discussed (and is 
obvious), there can be no low-energy 2D plasmon modes 
in intrinsic graphene. 

Putting /X = and taking the long-wavelength limit 
ypq/ksT — >• 0, we get from Eqs. ([T]) and ^ the follow- 
ing expression for the finite-temperature {ksT/vFq ^ 1) 
intrinsic graphene polarizability function (We mention 
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FIG. 2: (Color online). Plasmon damping rate of intrinsic 
MLG. (a) presents the plasmon damping rate in the high tem- 
perature regime (fcsT ^ hvpq)- The dashed and solid lines 
correspond to the analytical results (given in Eq. ((8)1 ) and 
the numerical results, respectively, (b) Numerical results for 
h'y/{kBT) as a function of hvFq/{kBT). (c) 7p/aj as a func- 
tion of hvpq/ {ksT). The dashed and solid lines correspond 
to K = 1 and K = 5, respectively. 
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FIG. 3: (Color online). Results for extrinsic MLG with 
q/kp ~ 0.01. (a) presents the plasmon dispersion in the low 
temperature regime. The dashed and solid lines correspond 
to the analytical results (given in Eq. and the numerical 

results, respectively, (b) Numerical results for uop/Ep as a 
function of T/Tp- (c) Numerical results for /Ep as a func- 
tion of T /Tf- (d) Numerical results for ')/uJp as a function of 
T/Tf- The dashed and solid lines correspond to k = 1 and 
K = 5, respectively. 



that the complete analytical n(g,a;) at T = for both 
intrinsic and extrinsic case, i.e. intraband and interband, 
can be found in Ref. [l|): 



n(g,w) 



21n2 



knT 



16 _ q2^2^ keT 
which reduces to (in the limit uj 3> qvp as q — > 0): 

21n2g2 i 
Ii{q,uj) w ^ksT 



16 keT 



(5) 



(6) 



Putting Eq. in Eq. ([T|), and solving for the complex 
frequency defining the intrinsic plasmon, we get: 



= {4:\n2)rshvFqkBT 
7 = g^fc^ V(ln2)rsfii'fg 



(7) 
(8) 



We have restored K in Eqs. ([7]) and ^ for the sake of 
clarity and usefulness (and = / nhvp is the graphene 
fine-structure constant). 

Equations ([7|) and ([8]) define the long- wavelength (and 
necessarily finite temperature) intrinsic MLG plasmon, 
which is determined by the variables temperature and 
wavevector (and not by any carrier density n as in all 
ordinary plasmon modes). We note that Wp,7 — >■ as 
g 0, but obeying different power laws: uip ~ ,Jq con- 
sistent with the 2D plasmon behavior and 7 ^ q'^/'^. We 



also note that Wp '-^ a/T and 7 ^ \/\/T whereas ujp 



3 /2 . 1 . . 

and 7 ~ Ts with ^ giving the dependence on 
the background dielectric constant. 

We note that cjp/^ = 16 nk bT / (Tre'^q) from Eqs. ^ 
and ([5]) , and therefore the long- wavelength intrinsic plas- 
mon is well-defined as long as 



q < 



16KfcsT 



(9) 



which defines the condition for the existence of a well- 
defined long-wavelength intrinsic MLG plasmon. Thus, 
the intrinsic plasmon remains well-defined at long- 
wavelength for arbitrarily low temperature as long as one 
is probing wavevectors shorter than the critical wavevec- 
tor qc defined by 



9c 



WKksT 



(10) 



For q < qc, the MLG intrinsic plasmon exists as a well- 
defined long wavelength collective mode, and for q > q^ 
it is overdamped (i.e. 7 > Wp). 

Before providing our full numerical results for the MLG 
intrinsic plasmon for arbitrary q and T, we briefly com- 
pare the analytical results for intrinsic and extrinsic MLG 
plasmons, which were earlier considered in Refs. [l[ and 
[lOj. The T = plasmon dispersion for extrinsic (i.e. 
doped) graphene is given in the long wavelength limit 
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FIG. 4: (Color online). Results for extrinsic MLG with 
q/kp = 0.5. (a) Numerical results for uj^/Ep as a function 
of T/Tp. (b) Numerical results for 'y/Ep as a function of 
T /Tp. (c) Numerical results for /uip as a function of T jTp. 
The dashed and solid lines correspond to k = 1 and k = 5, 
respectively. 
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FIG. 5: (Color online). Results for extrinsic MLG with 
qjkp — 1.0. (a) Numerical results for ujp/Ep as a function 
of T/Tp. (b) Numerical results for "y/Ep as a function of 
T /Tp. (c) Numerical results for 'y/u)p as a function of T/Tp. 
The dashed and solid lines correspond to ft = 1 and ft = 5, 
respectively. 



by: 



Lip = {2rshvpqEF 



,1/2 



1/2 



(11) 



where n is the carrier density (with a Fermi level Ep — 
Hvf^/t^ 7^ 0). It is easy to obtain the low-temperature 
analytical result for the MLG plasmon dispersion by us- 
ing the finitc-tcmpcrature expansion for the chemical 

2" 

1 



potential: fi{T) w Ef 
Ep/ks. We get for T ^Tp: 



for T ^ Tp 



UJp ~ ( 2rshvpqEp 



^ 6 \Tp 



y 2X. 



1/2 



(12) 



Thus, ujp{T) = iOp{T = 0) 



1 - — 

^ 12 



which 



T 

is a small correction to the T = result. We note 
that the Landau damping for extrinsic plasmons at long 
wavelengths and low temperatures is exponentially sup- 
pressed, going as e""^/"^^ . We remark here that the rea- 
son that the finite-T extrinsic plasmon is exponentially 
weak by Landau damped (as T — > 0) whereas the corre- 
sponding intrinsic plasmon has power law (7 ~ l/VT) di- 
vergent Landau damping as T — ?> is that, by definition, 
the intrinsic plasmon is always in the high-temperature 
regime for any temperatures since Ep ~ kgTp = 
for intrinsic plasmons. (Below we will discuss the high- 
temperature limit for the extrinsic plasmon.) 



An interesting exercise (alluded to in the Introduction 
of this paper) is to ask whether the long-wavelength in- 
trinsic plasmon dispersion (i.e. Eq. ([7])) can be obtained 
from the corresponding extrinsic plasmon dispersion (i.e. 
Eq. (jlip ) by reinterpreting the doping carrier density n 
in Eq. (jlip as the thermally excited carrier density n{T) 
for the intrinsic case. The thermally excited carrier den- 
sity n(T) for intrinsic graphene with the Fermi level at 
the Dirac point {Ep = 0) is easily calculated to be: 



n{T) = J D{E)f{E)dE 



(13) 



where D{E) = AE/{2Trvp) is the graphene density of 
states. Integrating over the Fermi distribution function 
f{E) at temperature T with /i = wc get: 



Inserting Eq. ([T51) for n in Eq. ([TT|) we get: 



(14) 




rs{hvpq){kBT) 



1/2 



(15) 



Eq. ([T5|) has the same parameter dependence 
^rshvpqksT as in the correct intrinsic plasmon 
dispersion given by Eq. ([7]) with the only difference is 
that the prefactor in Eq. (fTSj) is 2y7r2/6 ~ 2.6 versus 
the prefactor in Eq. (O is 4 In 2 sa 2.8. Thus, the plasma 
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FIG. 6: (Color online). Plasmon dispersion as a function of 
q/kp for extrinsic MLG. (a) and (b) present uip/Ep in the 
long wavelength limit. The dashed and solid lines correspond 
to the analytical results (given in Eq. (|12|) ') and the numerical 
results, respectively, (a) is for T = Q and (b) is for T/Tf ~ 
0.1. (c) and (d) present ajp/_EF for different temperatures, (c) 
is for K = 1 and (d) is for k — 5. 
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FIG. 7: (Color online). Plasmon damping rate as a function 
oiq/kp for extrinsic MLG. (a) and (c) '^/Ef as a function of 
q/kF for different temperatures, (b) and (d) present 7/ajp as 
a function of q/kp- (a) and (b) are for n = 1; (c) and (d) are 
for K = 5. Note that the legend applies to all sub-figures. 



frequency is given by 1.67-\/rs hypqiksT)^/'^ and by 
1.60V»'..fiuFg(fcBT)i/2 in Eq. 

We now consider the high-temperature limit of the ex- 
trinsic plasmon dispersion for gated or doped graphene 
taking T ^ Tp{= Ep/ks)- The asymptotic high- 
temperature expression for Ii{q, uj) in doped graphene is 
given by (again for qvp 'ti uj): 



, 2 In 2 fl^ 
n(g,w) « -^fcsT ( 1 



+ 



g2 keT 



128(ln2)3r4 

1 - 



48fc|T2 



(16) 



We emphasize that Eq. ([T6|) is valid for extrinsic graphene 
in the limit of T > Tf and qvp -^uj. Using Eq. to 
solve for the plasmon modes in Eq. ([T]), we get: 



7 = 



y41n2 rshvpqksT ( 1 
nhvpqrs 



Tp 



128(ln2)3r4 



(17) 



SVkBT 



A direct comparison between Eqs. ([T7)) . ([T5)) and Eqs. ([7|), 
([5]) show that the T/Tp — > oo limit of the extrinsic plas- 
mon dispersion and broadening indeed agree with the 
corresponding intrinsic plasmon results in the leading 
order, as indeed it must. We mention, however, that 
this agreement is only in the T/Tp — )■ oo limit (i.e. the 
Tp = limit of the extrinsic situation). Thus, there is a 



correction to the leading-order extrinsic plasmon disper- 
sion in Eq. (|17p going as 0(T^/T^) which, by definition, 
cannot exist in the intrinsic plasmon dispersion where 
the leading-order dispersion comes entirely as y/qT with 
no correction term in temperature. All higher-order tem- 
perature corrections to the intrinsic plasmon dispersion 
occur in higher-order terms in the wavevector q. 

We now present in Figs. [T] and [2] our calculated numer- 
ical results for the intrinsic graphene energy dispersion 
and level broadening for arbitrary values of wavevector 
q and temperature T. (Our analytical results presented 
above in Eqs. (O and ([5]) are necessarily restricted to 
the hvpq <^ ksT regime.) In presenting our results, 
we find that there are only two dimensionless (indepen- 
dent) variables that completely characterize the intrin- 
sic graphene plasmon properties: k, and Hvpq/kBT. In 
Figs. [1] and [2l we use two values of k = 1 (suspended 
graphene) and 5 (graphene on boron nitride (BN) sub- 
strate) as representative examples of strongly (k = 1, i.e. 
rg « 2.2) and weakly {k = 5, i.e. rg ~ 0.4) interacting 
systems to show our numerical results for the numeri- 
cally calculated plasmon energy hujp and level broaden- 
ing h'y in units of ksT as functions of the dimensionless 
2D wavevector hvpq/ (ksT). Our numerical results solve 
Eq. ([1]) to obtain the complex solution with the real part 
being the plasma frequency and the imaginary part the 
broadening. We obtain Il{q, lo) at arbitrary temperatures 
numerically in order to solve for the plasmon modes at 
arbitrary temperatures and wavevectors. In Fig. [1] the 
plasma frequency is shown as a function of wavevector, 
both for the small- g regime and over a large range of q. 
The small-g results serve to verify the accuracy of our 



6 




\ (b) g/A:F=0.1 
O.OlSf K = l 



0.00 0.02 O.0t 0.06 0.08 0.10 

o.om 



Om 0.04 0.06 0.0 

q/<.-f 




0.000 1 

0.00 0.02 OM 0.06 0.08 O.IO 

q/kf 



FIG. 8: (Color online). Results for extrinsic MLG with 
T/Tf = 10.0. (a) presents plasmon dispersion in the long 
wavelength limit. The dashed and solid lines correspond to 
the analytical results (given in Eq. (IT7|) and the numerical 
results, respectively, (b) Plasmon damping rate versus q/kp- 
The dashed and solid lines correspond to the analytical results 
(given in Eq. {18])) and the numerical results, respectively, (c) 
^/ijJp versus q/kp- The dashed and solid lines correspond to 
the analytical results (given in Eq. (|17l) ) and the numerical 
results, respectively. 




FIG. 9: (Color online). Results for extrinsic MLG. (a) and 
(c) present plasmon dispersion versus T/Tf- The dashed 
and solid lines correspond to the analytical results (given in 
Eq. p7[) ) and the numerical results, respectively, (b) and (d) 
present plasmon damping rate versus T/Tf- (a) and (b) are 
for q/kp ~ 0.1. (c) and (d) are for q/kp = 0.5. The dashed 
and solid lines correspond to the analytical results (given in 
Eq. (118^ ) and the numerical results, respectively. 



asymptotic analytic result given in Eq. ([7]) . In Fig. [2] we 
depict our calculated plasmon broadening as a function of 
wavevcctor again for small-g (Fig. dja)) and cxtended-g 
(Fig.[2jb)) regions whereas in Fig. [2jc) we depict the di- 
mensionless ratio as a function of the dimensionless 
variable hypq/ksT . 

The most notable, and perhaps somewhat unexpected, 
feature of our numerical results in Figs. [T] and [3] is that the 
intrinsic plasmon mode remains well-defined, i.e. uip > 7, 
for all values of Hvpq/kBT with a shallow maximum 
around Hvpq ^ ksT manifesting a surprising nonmono- 
tonic behavior for both values of k in Fig.[2t^c). For k = 1 
(rs = 2.2), suspended graphene, the maximum value of 
oopf-f reaches 0.4, but for k = 5 (r^ = 0.4), graphene 
on BN, the maximum value of ujp/^ is below 0.1. Thus, 
the intrinsic MLG plasmon should, in principle, be ex- 
perimentally observable, particularly on substrates with 
large dielectric constant where 7/0;^ <C 1 making the 
Landau damping problem fairly irrelevant. Our results 
in Figs. [1] and [2] also indicate that the leading order for- 
mula of Eqs. ^ and ([8]) remain reasonably well- valid for 
arbitrary values of hvpq/kBT . 

Although our focus in the current work is the intrinsic 
Dirac point plasmon mode for undoped graphene, it is 
useful to compare the temperature dependence of the ex- 
trinsic plasmon in doped graphene with that of intrinsic 
graphene, particularly since the temperature dependence 
of extrinsic graphene plasmon has not much been studied 



in the literature. Understanding temperature dependent 
plasmon dispersion and damping of doped graphene is 
also relevant here since extrinsic and intrinsic graphene 
plasmons become the same at very high temperatures 
(T ^ Tp). We therefore provide a large set of finite- 
temperature results for extrinsic plasmon dispersion and 
damping, comparing them with our derived analytical 
low-and high-temperature results and with the corre- 
sponding intrinsic plasmon results. Our motivation for 
such a detailed finite-temperature RPA study of extrin- 
sic graphene plasmons comes partially from the fact that 
temperature could in principle be used (in addition to 
wavevector and/or carrier density) to tunc the plasmon 
energy in graphene (particularly at lower carrier densities 
and higher temperatures where T /Tp is not necessarily 
extremely small) , a fact which has not been much appre- 
ciated in the literature. 

In showing our full numerical solutions for u}p{q) 
and 7(g) using Eq. ([T]) [and the full finite-temperature 
n(g, uj)] for extrinsic graphene, the first problem we face 
is that there are far too many independent variables 
(i.e. q, T, n, k) which determine the plasmon properties. 
Since three of these variables are independent continu- 
ously tunable experimental variables (i.e. q,T,n), we 
need four-dimensional plots (for several values of k, i.e. 
for different substrates) or perhaps even five-dimensional 
plots showing ujp and 7 as functions of q, T, n, and k. A 
significant simplification arises from using kp(~ \/Tm) 
and Ep(= hvpkp — hvp^/wn) as the unit of wavevec- 
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FIG. 10: (Color online). Results for extrinsic MLG with 
q/kp = 0.5 over a wide T /Tp = — 10 range, (a) Numerical 
results for ujp/ Ep as a function of T /Tp. (b) Numerical re- 
sults for 'y/Ep as a function of T /Tp. (c) Numerical results 
for "i/ujp as a function of T/Tp. The dashed and solid lines 
correspond to k = 1 and k = 5, respectively. 
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FIG. 11: (Color online). Results for extrinsic MLG with 
q/kp = 1.0 over a wide T /Tp = — 10 range, (a) Numerical 
results for uip/ Ep as a function of T /Tp. (b) Numerical re- 
sults for 'y/Ep as a function of T /Tp. (c) Numerical results 
for as a function of T /Tp. The dashed and solid lines 

correspond to «: = 1 and «: = 5, respectively. 



tor and energy, respectively so that the carrier density 
shows up implicitly as a scaling variable rather than ex- 
plicitly, eliminating one variable. We also show results 
only for two values of the background dielectric constant 
K = 1 (suspended graphene) and 5 (graphene on h-BN 
substrates) corresponding to = e^/{hvpK) = 2.2 and 
0.4, respectively. Thus, we present all our numerical re- 
sults for ujp/ Ep and ^ /Ep as functions oiq/kp and T /Tp 
for /t = 1 and 5 in the following. A particular goal of the 
presented numerical results for arbitrary q and T is com- 
parison with our analytical low- and high-temperature 
results in Eqs. ((T^ and ((T7)) / ([T5|) . respectively. Since we 
will be presenting a very large number of figures, we do 
not discuss all the figures individually in the text below, 
instead only highlighting the important salient features. 
We provide detailed figure captions in the figures them- 
selves which should be self-explanatory. 

In Figs. [5][S] we present results function of T /Tp 
for fixed q/kp = 0.01;0.5;1. The extrinsic plasmon en- 
ergy for fixed q generically shows a non-monotonic depen- 
dence on temperature with a shallow minimum around 
T ~ Q.ATp. This arises from the fact that the degenerate 
system (T <C Tp) has a plasma frequency decreasing with 
increasing temperature (according to Eq. (|12p ) whereas 
the nondegenerate system has the plasma frequency in- 
creasing (as ^ T for T ^ Tp) with temperature. The 
broadening 7 is suppressed exponentially for small T /Tp 
except for large q{^ kp) where intraband Landau damp- 
ing starts playing a role, particularly for larger values 
of Ts (smaller n). The analytic formula, Eq. for 
the low-temperature plasmon dispersion (and damping) 



seems to work very well (somewhat surprisingly) all the 
way to T/Tp ^ 0.4, i.e. all the way to the shallow 
minimum in u)p{T). According to Figs. the plas- 
mon damping manifests a shallow maximum around the 
same value of T/Tp ('^ 0.4) where the plasmon energy 
is a maximum, and thus "f/ujp shows a generic peak for 
T/Tp ^ 0.4 (which is much sharper for larger values). 
In general, we find j/tUp < 1 for smaller q values in the 
T/Tp < 1 regime — thus the plasmon is well-defined for 
T < Tp. 

In Figs. [BJH wc show the plasmon energy and broad- 
ening at various finite T values as a function of wavevec- 
tor — thus, Figs. [B][5] are effectively complementary to 
Figs. [3][5l First, we note that the analytic results are 
essentially in exact agreement with the full numerical re- 
sults upto q ^ 0.5kp. In Fig. [71 we show that 7 < is 
well-satisfied for all q values well up to T ~ Tp. In gen- 
eral, the broadening is exponentially suppressed at low q 
and low T, but increasing either qovT eventually leads to 
intraband and interband Landau damping. In Fig. [71 the 
onset of the intraband Landau broadening, where the ex- 
trinsic plasmon dispersion enters the intraband electron- 
hole single particle excitation continuum even at T = 0, 
can clearly be seen around q ^ 0.4 — 0.7kp, whereas for 
higher T values, the extrinsic plasmon can decay even at 
long wavelength due to interband electron-hole excitation 
process. 

Whereas in Figs. [3|[7] we focus on low-temperature 
(T < Tp) extrinsic plasmon dispersion and damping, we 
now present in Figs. MT^ higher-T(> Tp) results. The 
higher-temperature (T > Tp) extrinsic plasmon results 
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FIG. 12: (Color online). Results for extrinsic MLG with k = 
5. (a) Numerical results for ujp/Ep as a function of q/kp- 
(b) Numerical results for "j/Ep as a function of q/kp. (c) 
Numerical results for 'y/up as a function of q/kp. The solid, 
dashed, dotted and dot-dashed lines correspond to T/Tp = 
0, 1, 5, and 10, respectively. 
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FIG. 13: (Color online). Results for extrinsic MLG with k = 
1. (a) Numerical results for Up/Ep as a function of q/kp. 
(b) Numerical results for 'y/Ep as a function of q/kp. (c) 
Numerical results for as a function of q/kp. The solid, 

dashed, dotted and dot-dashed lines correspond to T /Tp = 
0, 1, 5, and 10, respectively. 



are relevant for understanding intrinsic plasnion behav- 
ior in graphene since, as emphasized in our analytical 
theory (see Eqs. (|17p. ((TH)) . and the discussions following 
their derivations), the leading-order (in Tp/T) results for 
both plasnion energy and damping for intrinsic and ex- 
trinsic plasmons are the same for T ^ Tp. Physically, 
the reason for this is obvious: For T ^ Tp, the thermal 
inter-band electron-hole excitations dominate the collec- 
tive behavior over the contribution by the doped carriers 
even for extrinsic doped graphene (of course, the actual 
temperature scale needed to satisfy the T ^ Tp condi- 
tion increases as \/n with increasing the doping density) . 

One of the most important as well as interesting as- 
pects of the numerical results shown in Figs. [8llT3lis the 
great quantitative accuracy of our leading order analytic 
high-temperature (T ^ Tp) results (i.e. Eas. p7|) and 
(|18p ) for extrinsic plasnion dispersion and damping as 
compared with the full RPA finite-temperature numeri- 
cal solutions for Wp(g) and 7(T). In particular, our ana- 
lytical theory seems to hold very well all the way down 
to T > Tp although the analytic theory represents an 
asymptotic expansion in Tp/T. This remarkable reliabil- 
ity of the high-temperature analytic theory (for T ^ Tp) 
as well as that of the low-temperature {T ^ Tp) analytic 
theory (derived as an asymptotic expansion in T/Tp), 
which we discussed in the context of Figs. |3]l7] above, 
implies that the analytical finite-temperature theory de- 
veloped in the current work could be extremely useful 
for experimental works in graphene plasmonics with no 
need for the full numerical solution of the RPA theory 
which is, in fact, quite complex and demanding at fi- 



nite temperatures since the finite-temperature graphene 
polarizability Il{q,iL);T) does not have any simple ana- 
lytical form and must be carefully calculated through a 
numerical integration at each value of T. If necessary, 
one could easily develop a numerical interpolation scheme 
between our low-temperature and high-temperature an- 
alytical theories (e.g. using a suitable Pade approxima- 
tion) which should provide a reasonable and quantita- 
tively accurate theory at arbitrary temperatures. Since 
the high-temperature analytical extrinsic plasmon the- 
ory (i.e. Eqs. PT)) and dH])) essentially agree with the 
intrinsic plasmon results, the T 3> Tp results provided 
in Figs. [51fT51 could be construed as numerical results for 
the intrinsic plasmon as well. 



B. Double layer graphene 

Collective modes of two (i.e. "double") parallel 2D 
graphene layers (along the x—y plane) separated by a dis- 
tance "d" in the third direction (z— direction) were first 
theoretically considered by Hwang and Das Sarma^S and 
later by other author a^^i^i" — . In the current work, we fo- 
cus on the double-layer system considering their plasmon 
modes at finite temperatures, when one (or both of the 
layers) is (are) intrinsic or undopcd. Thus, our work is 
the finite-temperature generalization of the Hwang-Das 
Sarma work, concentrating on intrinsic plasmons in un- 
doped double layers. 

The collective modes of a double-layer system is 
obtained by diagonalizing the 2x2 determinantal 
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equation's defined in Eq. which gives: 

(1 - yiiiinii)(i - y2222n22) - v^\,^nnn22 = o (i9) 

where 11;; (with / = 1,2) is the polarizabiUty of the Itli 
layer, Vuu is the Coulomb interaction (i.e. 2Tre'^ / {kiq)) 
in the Ith layer and Vuiw is the Coulomb interaction be- 
tween electrons in the / and I' layers. We ignore here 
(and throughout this paper) the possibility of electron 
hopping between the two layers which is an excellent ap- 
proximation for graphene. For our double-layer system 
we have: 



Vim — V2: 



and 



Vl2U 



2^ 
Kq 



Kq 



(20) 



(21) 



where we have assumed both layers to be submerged 
in the same background dielectric with k as the 
common background lattice dielectric constant. (A 
generalizatioii^ to the situation with ki, K2 is straight- 
forward, but will not be considered in our work since it 
will add more parameters to a problem which already has 
far too many variables.) 

Using the known analytical expressions for 11 = ITn = 
1122 in the long-wavelength limit for intrinsic graphene 
(Eqs. (O and (H])), we can solve Eq. to get the follow- 
ing two coupled long-wavelength collective modes {ujop 



and 



for the double-layer system when both layers 



are intrinsic graphene: 



lop 
lap 



^{%\n2){hvFq){kBT)rs 
^{4hi2){rsdq){hvFq)kBT 

' '^V(ln2)r shvpq 



irrshvpq^d 



\/ (In 2)rshvFq'^d 



(22) 
(23) 

(24) 
(25) 



Equations (|22p - (P5)) define the long wavelength collective 
modes and their damping for an intrinsic double-layer 
graphene system where both layers are undopcd (and the 
Fermi level in both layers is sitting at the Dirac point). 

Wc first discuss the implications of our derived 
(Eqs. ([22]) - ([25]) ) analytical results for double-layer 
graphene intrinsic plasmons. It sounds crazy that two 
undoped graphene layers (i.e. no free carriers whatso- 
ever) can have two low-energy collective modes {ujop and 
ujap above) when they are proximate to each other, but 
it is apparently true. One of these modes, ujop, is noth- 
ing other than the combined collective intrinsic plasmon 
mode of each independent layer with = Sw^, where 
Up is the intrinsic plasmon frequency of a single graphene 
layer as given in Eq. ([7]). Thus, Wop is simply the in-phase 



intrinsic plasma oscillation of the two intrinsic plasmons 
in the two layers. The uJop mode is sometimes referred 
to as the "optical plasmon" mode^S of the 2— component 
(i.e. two layers) double-layer system since it involves the 
in-phase collective charge density oscillation of the two 
layers (analogous to an optical phonon mode in a lattice). 
The second mode, uJap, which has no analog in the single 
layer system, is the acoustic plasmon modoSS where the 
charge density oscillates out of phase between the two 
layers. 

The aJop(oc ^/q) obviously has the same dispersion 
as the single layer plasmon whereas the ujap (x q has 
an acoustic dispersion linear in q at long wavelengths. 
Both modes have the basic intrinsic plasmon property of 
LOp oc Vt as expected; and vanishes in the T — )■ limit. 
The broadening (Eqs. ([^ and (^5])) has a higher-order 
g— dependence, q^^^{q^) for Wopi^ap), respectively, ensur- 
ing that both modes are well-defined collective modes in 
the long-wavelength limit. In particular, we have: 



^ op hop = 



■ne^q 

ap lap — 2 2 J 



(26) 
(27) 



Thus, the optical plasmon mode of the double layer sys- 
tem has very similar behavior for ajp/7 as in the corre- 
sponding single layer intrinsic plasmon case except that 
the ratio of Wp/^ is a factor of 2 smaller for the dou- 
ble layer case than the single-layer case. The Wop mode 
therefore remains well-defined (i.e. > 7op) down to a 
cut-off wavevcctor 



qoc = 



ne 



(28) 



with ojop > Jop for all q < qoc- 

The situation for the intrinsic acoustic plasmon ujap 
is, however, qualitatively different since uJap oc q (rather 
than y/g as for ujop) and thus it decreases fast as g — >■ 0. 
The cut-off wavevector for a well-defined Uap mode is: 



KksT 
dne^ 



(29) 



> 7ap for all q < qa 



with I 

Wc point out that the analytical results for the acous- 
tic plasmon given above (Eqs. (|23p and (|25t ) apply only 
when ksT > hvp/ (4 In 2 rgd) in order to satisfy the crite- 
rion uj > hvpq used in the expansion of Il{q,uj) to derive 
the long-wavelength plasmon dispersion relation. 

In Figs. [Ml and [TSl wc provide our full RPA numerical 
results for the double-layer intrinsic plasmon modes and 
compare them with our analytical results (Eqs. (^^ - ((25)) ) 
obtained above. To keep the number of presented figures 
tractable we only change q and T for showing our results. 
In general, the theory no longer scales with hvpq/kgT 
as the corresponding single-layer intrinsic problem docs 
because of the presence of the layer separation d in the 
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problem. For a fixed "d" and "T" , however, we can still 
show results using fc^T as the energy unit remembering 
that these double-layer results of Figs. [M] and [15] apply 
only for fixed results of k, d, T as shown in the figure (but 
for varying q). 

In Fig.[T4]we show the coupled intrinsic plasmon modes 
of double-layer graphene for small values of hvpg/ksT 
where our analytical expressions derived in Eqs. (|22p - 
(|25|) are essentially exact. The uJop (u^ap) modes show 
the expected y/q (q) dispersion, and typically ujop > lop 
{ujap > lap)- 111 Fig- US] we depict the typical plasmon 
dispersion and damping for the coupled modes for three 
values of K, d, T over a broader range of hypq/ksT, find- 
ing that for layer values of d{~ 300 A) and smaller k(= 1) 
the damping could be quite large. 
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FIG. 14: (Color online). Results for an intrinsic double-layer 
graphene system in the small-g regime with k = 1, d = 300 A 
and T = 100 K. (a) Plasmon dispersion versus hvpq/ {ksT). 
(b) Plasmon damping rate versus Tmpq/ {ksT). The dashed 
and solid lines correspond to the analytical results (given in 
Eqs. (|22|l - (|25|l ') and the numerical results, respectively. 



Next, we consider the double-layer graphene plasmons 
when one layer is doped ("extrinsic") and one undoped 
("intrinsic"), i.e. one layer has carrier density n 7^ and 
the other has n = at T = 0. The long wavelength and 
low temperature RPA collective modes of such a mixed 
intrinsic-extrinsic graphene double-layer system are eas- 
ily derived to be given by: 



2rsh^VpqkF ( 1 + (2 In 2) 



T 
7> 



7r2 t2 

6 r| 



1/2 



(30) 



lop 



lap 



I (21n2)T 4(ln2)2T2n ^ 
na2)r,hvpq^dkBT\l- ^ ' + (31) 



n^r^H^Vpq'^Ep 
128fc|T2 



1 + (2 In 2) 



Tp 
T 

3> ' 



rp2 
^ P 



1 1/2 



(32) 



t2 In 2r'ih?v\.q^d? 



SkgT 



[l-(101n2)^+60(ln2)2g^3) 



T2- 



Here = e'^/{Khvp) refers to both layers, and kp = 
{imy^^ and Tp = Ep/ks = hvpkp/kB refer to the 
doped extrinsic layer. 

While the above results are valid in the T/Tp <^ 1 
(as well as leading order in q) limit, we can also ob- 
tain the high-temperature {T/Tp 3> 1) asymptotic an- 
alytical results to be exactly the same as those given in 
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FIG. 15: (Color online). Numerical results for an intrinsic 
double-layer graphene system, (a) and (b) Plasmon disper- 
sion versus faipqj (ksT); (c) and (d) Plasmon damping rate 
versus hvpq/{kBT); (e) and (f) l/ujp versus hvpq/{kBT). (a), 
(c) and (e) are for the optical plasmon mode, (b), (d) and (f) 
are for the acoustic plasmon mode. For «; = 5, d = 40 A, and 
T = 20 K, the acoustic plasmon mode ujap is degenerate with 
the boundary of the single particle excitation region. Note 
that the legend applies to all sub-figures. 



Eqs. (|^ -(P5 |) for double- layer intrinsic graphene. This 
is, of course, expected since in the T ^ Tp limit, there is 
no difference in the leading order between intrinsic and 
extrinsic graphene. 

The mixed double-layer intrinsic-extrinsic graphene 
system depends in a complicated manner on a large num- 
ber of independent parameters: q,n,T,d, k. The plas- 
mon modes now really depend on all five of these param- 
eters (plus the value of vp which in principle is also a 
free parameter of the theory). We refrain from overload- 
ing the readers with a large number of results varying all 
five parameters freely. Wc provide some representative 
plasmon dispersion and damping results in Figs. I16II20I 
showing the full numerical RPA solutions for the plas- 
mon energy and damping for the mixed double-layer sys- 
tem, emphasizing that our low-temperature (high tem- 
perature) analytical results seem to work very well in the 
T < Tp (T > Tp) regimes, providing an easy and effec- 
tive way for quick calculations of the plasmon dispersion 
and damping in double-layer intrinsic-extrinsic plasmon 
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systems. The captions in each figure (in Figs. [T6ll20)) 
clearly describe in details the parameter values and the 
numerical results for the mixed intrinsic-extrinsic double 
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FIG. 16: (Color online). Results for a mixed double-layer 
intrinsic-extrinsic graphene system with ni — IQ^^ cm^^, 
'n,2l'n\ = 0.0, d — 300 A and k — \. (a) ^jp/Ep versus q/kp 
for T/Tf = 0.1. m-y/Ep versus q/kp for T /Tf = 0.1. (c) 
uJj,/Ef versus T /Tf for q/kp = 0.01. [A)^ / Ef versus T /Tf 
for q/kp = 0.01. The dashed and solid lines correspond to the 
analytical results (given in Eqs. (|30 p -(|33 [l l and the numerical 
results, respectively. 



Finally, wc conclude this subsection with a brief discus- 
sion of the finite-temperature double-layer extrinsic plas- 
mon system when both layers are doped. Since the T = 
case for this system was considered in details by Hwang 
and Das Sarma2£ with some follow-up T ^ Q double- 
layer plasmon calculations in the literatur o^^'"^^" — , we 
only provide some brief analytical results for the sake 
of completeness (and comparison with our intrinsic plas- 
mon double layer results). Wc mention that the extrinsic 
double-layer plasmons arc in principle determined by six 
independent parameters (rii, 712, T, q, d, n), and providing 
complete numerical results here will simply make our pa- 
per far too long. 

Using the asymptotic forms for the graphene polariz- 
ability in the long wavelength limit, we get the following 
analytical formula for the low-temperature (T <?C Tp) 
plasmon dispersion and damping in the leading order: 



l2r,h^vj.q{kFi+kF2) 1- 
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FIG. 17: (Color online). Numerical results for a mixed 
double-layer intrinsic-extrinsic graphene system over a wide 
q/kp =0 — 1 range with n\ = 10^^ cm^^, n^/nx = 0.0, 
d = 300 A and k = 1. (a) ujp/Ep versus q/kp for different 
temperatures, (b) "y/Ep versus q/kp for different tempera- 
tures, (c) ')/uip versus q/kp for different temperatures. For 
T /Tp = 0, the acoustic plasmon mode Uap is degenerate with 
the boundary of the single particle excitation region. Note 
that the legend applies to all sub-figures. 



(35) 



where fcFi,2 = ksTpia = Efi,2 = fivpkFia 

depend on the doping carrier density ?7i.2 in the two lay- 
ers, and T <C Ifi_2 is assumed in obtaining Eqs. ((34)) 
and (j35p . The corresponding low-temperature long- 
wavelength damping 7 for both plasmon modes is expo- 
nentially suppressed as in the corresponding single-layer 
extrinsic plasmon case. 

We can also carry out the high-temperature T ^ Tpi 2 
asymptotic expansion of 1X1.2 to derive the corresponding 
high-temperature results and wc get precisely Eqs. (|22p - 
(p5)) in the leading order in TF1.2/T, i.e. Tp drops out 
in the leading order, leaving us precisely the intrinsic 
double-layer plasmon results, as expected. The precise 
agreement between the extrinsic and the intrinsic results 
arises in the leading order in Tp/T because the finite- 
T chemical potential in the extrinsic case goes as /i « 

We just show three representative sets of numerical 
results for double-layer extrinsic graphene in Figs. [2Tll23l 
In Fig. (211 the small-(j and small-T numerical results are 
shown manifesting excellent agreement with our T <^ Tp 
analytical results. In Fig. [22] we show the results for 
the same parameters in a much more expanded scale of 
q/kpi. In Fig. [^niwc show the temperature dependence 
at fixed q. 
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FIG. 18: (Color online). Numerical results for a mixed 
double-layer intrinsic-extrinsic graphene system over a wide 
T/Tf = 0-1 range with ni = 10^^ cm~^, n2/ni = 0.0, 
d = 300 A and k = 1. (a) ujp/Ep versus T /Tf for different 
q/kp- (b) 'y/Ep versus T/Tf for different q/kF- (c) 
versus T /Tf for different q/kp- Note that the legend applies 
to all sub-figures. 
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FIG. 19: (Color online). Numerical results for a mixed 
double-layer intrinsic-extrinsic graphene system over a wide 
q/kp =0 — 1 range with ni = 10^^ cm~'^, n2/ni = 0.0, 
d = 40 A and «; = 5. (a) ujp/Ef versus q/kp for different 
temperatures, (b) "y/EF versus q/kp for different tempera- 
tures, (c) ^/uip versus q/kp for different temperatures. For 
T/Tp = 0, the acoustic plasmon mode uiap is degenerate with 
the boundary of the single particle excitation region. Note 
that the legend applies to all sub-figures. 



C. Bilayer graphene 

Finally, we very briefly consider the theory for intrinsic 
plasmons in undoped bilayer graphene (BLG) for the sake 
of completeness. For our purpose, we assume the bilayer 
graphene to have parabolic chiral band structur o^^'^^" — 
(characterized by an cfFcctive mass m rather than a ve- 
locity Wi?) in contrast to the linear chiral band structure 
of MLG (which is what we have considered so far in this 
work) . 

The BLG 2D polarizability function is given by^: 



--y — 



,/s,k fs'M' 



IT] 



(36) 

where k' = k + q. s,s' = ±1, e^^k = sh^'k^/{2m) and 
Fs.s' = (1 + cos2^^)/2 is the chiral wavefunction overlap. 
The functions "/" in Eq. ([55)) are the Fermi distribution 
functions at temperature T. The bilayer dynamical po- 
larizability first calculated by Sensarma et al^ who also 
first obtained the bilayer graphene plasmon dispersion 
for the extrinsic (doped) system. Blow we give the ana- 
lytical results for intrinsic BLG plasmons in the absence 
of any doping. 

In the long wavelength limit {q — > 0) and at high 
temperatures h^q^/{2m) <C ksT, we have the follow- 
ing asymptotic formula for BLG n(g, w) in the intrin- 
sic {n = 0) regime (from a direct expansion of Eq. (|36p 



above): 



^ ^ 4 In 2 ^2 q 
Il{q,uj) « -ksT + i 



(37) 



Using Eq. ([T]) for the plasmon oscillation, we get the mode 
dispersion and damping for intrinsic BLG to be: 



8e2(hi2)(fcBT)g 



7 = 



ne q 



v/(ln2)e2g 



2Ky/2KkBT 

This gives for the BLG intrinsic plasmon: 
, SkBTn 

^vh = — — 



(38) 
(39) 

(40) 



For well-defined modes, we must have a;p/7 > 1, leading 
to the condition that: 



kBT> 



2 

ire q 
8k 



(41) 



Thus, similar to the MLG intrinsic plasmon although 7 oc 
1/\/T, there is a well-defined intrinsic plasma mode at 
long wavelength for arbitrarily low temperatures. 

Comparing the BLG intrinsic plasmon dispersion given 
in Eq. psp with the corresponding MLG expression given 
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FIG. 20: (Color online). Numerical results for a mixed 
double-layer intrinsic-extrinsic graphene system over a wide 
TjTp =0-1 range with m = 10^^ cm"^, n2/ni = 0.0, 
d = 40 A and /t = 5. (a) LUp/Ep versus T/Tp for different 
q/kp- (b) "i/Ep versus T /Tp for different q/kp. (c) "f/ujp 
versus T /Tp for different q/kp. Note that the legend applies 
to all sub-figures. 



FIG. 21: (Color online). Results for a double-layer extrinsic 
graphene system with ni = 10^^ cm~^, n2/n\ = 0.25, d = 300 
A. (a) and (b) u)p/Ep\ versus q/kp\. (c) and (d) ojp/Ep\ 
versus T /Tp\. The dashed and solid lines correspond to the 
analytical results (given in Eqs. (|34|) - H35I )) and the numerical 
results, respectively. 



in Eq. ([7]), wc sec that the MLG and the BLG results 
are identical except for an extra factor of 2 (8 In 2 in 
Eq. dSHl) and 4 In 2 in Eq. ^) in the BLG case. (Re- 
member that Ts = /(Khvp) making Eq. ([7]) equivalent 
to ujp ~ ^ (4 In 2)e'^kBTq/ K.) Thus the intrinsic plasmon 
frequency, in sharp contrast to the extrinsic plasmon fre- 
quency, is independent of the Fermi velocity (MLG) or 
the effective mass (BLG) and depends only on a single 
material constant k, the background dielectric constant. 
Wc note that the expressions for the intrinsic plasmon 
damping is also very similar for MLG (Eq. ([5])) and BLG 
(Eq. (PP]) ) — again except for a factor of VS, the two ex- 
pressions are identical (and independent of vp or m). 



the double-layer BLG system: 

ojop = \/(161n2)gfcsT— (42) 

V K 

oJap = \ {8ln2)q^dkBT— (43) 

V K 

These formula for the intrinsic double-layer BLG plas- 
mons are essentially identical to the corresponding 
double-layer MLG intrinsic plasmons (Eas. (P^ -(P5 |) ) ex- 
cept for numerical factors. Again, these plasmon mode 
dispersion and damping do not depend on the effective 
mass or the Fermi velocity and are universal properties 
of the intrinsic system. 

III. DISCUSSION 



We can also obtain the analytical results for the in- 
trinsic plasmons in the double-layer BLG system com- 
posed of two BLG layers separated by a distance "c?" . 
Solving the 2x2 determinantal equation (Eq. ([T5|) ) for 
the double-layer BLG system, we get the following an- 
alytic leading-order results for the intrinsic plasmons of 



Given the very large number of systems considered in 
this work (MLG, BLG, intrinsic, extrinsic, single layer, 
double-layer, • • • ) and the numerous analytical (as well as 
numerical) results presented in Sec|TTl it is useful to sum- 
marize all the analytical results in Tables HI and HIl so that 
the similarities/difFerenccs/connections among the vari- 
ous derived results for plasmon modes (and their damp- 
ing) become apparent. Some of the results in Tables |T] 
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FIG. 22: (Color online). Results for a double-layer extrinsic 
graphene system with ni = lO^'^ cm~^, n2/n\ = 0.25, d = 300 
A and hi = l. (a) cOp/Epi versus q/kpi for different tempera- 
tures, (b) -y/Epi versus q/kpi for different temperatures, (c) 
') /ujp versus q/kpi for different temperatures. Note that the 
legend applies to all sub-figures. 



FIG. 23: (Color online). Results for a double-layer extrinsic 
graphene system with ni = 10^^ cm~^, 712/711 = 0.25, d = 300 
A and K = 1. (a) uip/Epi versus T /Tpi for different q/kpi- 
(b) '/jEpx versus T jTpx for different q/kpi. (c) versus 
T /Tpi for different q/kpi- Note that the legend applies to all 
sub-figures. 



and [ni were obtained in the literature before, but our 
emphasis in this work is on intrinsic plasmon modes and 
their temperature dependence. Results in the last row of 
Tables HI and HIl are for the ordinary (non-chiral) 2D elec- 
tron gas systems (as occurring, for example, in semicon- 
ductor quantum wells) , which are by definition extrinsic 
systems since the large band gap (between conduction 
and valence bands) in the semiconductor ensures that 
for 71 = 0, there is no plasmon mode in the system. 

The key issue to be discussed in this section is the fea- 
sibility of an experimental observation of the graphene 
intrinsic Dirac point plasmon, which will be the mani- 
festation of a qualitatively new phenomenon since such 
intrinsic plasmons at zero doping density in a charge neu- 
tral system is essentially impossible in any 2D (or 3D) 
semiconductors where the existence of free carriers neces- 
sitates doping by some means. (As an aside we mention 
that even in a very narrow gap semiconductor, e.g. In 
As or InSb. the band gap is > 100 meV, which implies 
that T > 2000 K is necessary for any appreciable thermal 
population of free carriers making it impossible to study 
any kind of "intrinsic plasmons" in undoped semiconduc- 
tors.) 

There has been a great deal of recent experimental in- 
terest in studying the properties of graphene Dirac point 
(i.e. the charge neutrality point) in both MLG and BLG 
systems^"—. Various exotic quantum phases^^Ti^i and 
non-Fermi liquid behavior— are theoretically predicted 
at the graphene Dirac point, and thus the possibility 
of understanding the Dirac point behavior through the 
observation of intrinsic plasmon properties is both in- 



teresting and intriguing. For example, our theory of 
intrinsic Dirac point plasmon modes developed in this 
work is based entirely on the random phase approxi- 
mation assuming a generic Fermi liquid ground state, 
and consequently, any observed qualitative departure in 
the experimentally observed Dirac point plasmon behav- 
ior (e.g. a completely different temperature dependence 
compared with that given in Tables. |I]and|lT| would imply 
a failure of RPA indicating the fundamental and quali- 
tative importance of interaction effects or perhaps even 
the emergence of a new spontaneously symmetry-broken 
ground state as has been predicted theoretically in the 
Iitcrature22r-i2^. Since collective modes (e.g. Goldstone 
modes, Higgs bosons, zero sound modes, etc.) typically 
tell us a great deal about the fundamental nature of the 
field theoretic vacuum (i.e. the ground state of the sys- 
tem) , the study of graphene intrinsic plasmons could turn 
out to be a very useful route to understanding the nature 
of the Dirac point. A good example of the possible use- 
fulness of intrinsic plasmons could be in the determina- 
tion of whether the nonperturbative aspects of electron- 
electron interaction induce a chiral anomaly in graphene 
leading to the spontaneous formation of an energy gap 
at the Dirac point as has been predicted in some lat- 
tice quantum Monte Carlo simulationiS^iiS^iiSl. Such an 
energy gap at the Dirac point, if it exists in the experi- 
mental samples, would show up as an exponential (going 
as e~^l^ where A is the induced gap) suppression of the 
intrinsic plasmon energy (i.e. the simple power laws in 
temperature derived in our RPA theory will fail qualita- 
tively), which could be directly experimentally observed 
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thus validating (or not) the existence of a Dirac point 
chiral anomaly. 

Thus, an experimental investigation of intrinsic Dirac 
plasmons is highly desirable. It may appear to be a 
straightforward, even a trivial, task to study intrinsic 
Dirac plasmons experimentally by carrying out plas- 
mon experiments keeping the Fermi level fixed at the 
charge neutrality point by appropriately tuning the gate 
voltage in gated graphene. Transport measurements in 
gated samples routinely enable sitting precisely at the 
charge neutrality point since the graphene conductiv- 
ity has a typical "V" or "U" shape as a function of 
gate voltage with the conductivity (resistivity) minimum 
(maximum) being located at the nominal Dirac point. 
The serious complication is, however, disorder since it 
is well-knowni22r2ii that the charge neutrality point as 
determined by transport or other experiments is not 
the theoretical Dirac point because of the existence of 
electron-hole puddles in the system. In particular, ran- 
dom charged impurities in the graphene environment be- 
come qualitatively important in the n < Ui regime (where 
Hi is the impurity density and n the carrier density) , and 
drive the system into a highly inhomogeneous state with 
randomly distributed puddles of electrons and holes dom- 
inating the landscape. This puddle-dominated highly in- 
homogeneous (i.e. carrier density fluctuates spatially) 
regime has been show n^°^i^°^i^^^i^^^ to correspond to the 
so-called minimum conductivity plateau in the graphene 
transport data (i.e. the bottom of the "V" or "U" in the 
conductivity versus gate voltage plot). Thus, the charge 
neutrality point in transport corresponds only to the to- 
tal charge density in the whole sample being zero, and 
not to the absence of electrons and holes in the system. 
The puddle-dominated regime around the charge neutral- 
ity point should be best thought of as a random spatial 
variation in the Dirac point with respect to the spatially 
constant chemical potential (as controlled by the exter- 
nal gate voltage) where at each point in space the Dirac 
point is either above ( "hole regions" ) or below ( "electron 
regions") the chemical potential. Within the mean field 
theoryi22ii£^iii^"-ii^ the system has a finite fiuctuation 
induced carrier density {n* ~ n^) even at the Dirac point 
(i.e. the charge neutrality point) because of the Coulomb 
disorder induced electron- hole puddles. 

The existence of electron-hole puddles means that the 
Dirac point is ill-defined experimentally upto a carrier 
density of n* (with n* being determined by the details 
of the random charged impurity configuration in the sys- 
tem, but typically n* ^ rii), and n* can be approximately 
estimated experimentally by looking at the size of the 
conductivity minimum plateau regiortii^. The lack of the 
precise existence of the Dirac point has obvious implica- 
tions for the intrinsic plasmon which we discuss below. 

Electron-hole puddles make it impossible, as a matter 
of practice, to explore the precise Dirac point in graphene 
(or other similar materials), and therefore, our assump- 
tion of n = at the Dirac point is no longer applicable 
in our consideration of intrinsic plasmons. This, how- 



ever, should not be construed as a complete disaster for 
the observation of the intrinsic plasmon since the Dirac 
point, being a set of measure zero, would naturally be 
difficult to approach in any experiment^ii^, and all exper- 
imental claims of studying Dirac point phenomena are 
suspect because of the fragile and unstable nature of this 
"measure-zero" fixed point. Any experimental technique 
requiring the precise placement of the Fermi level at the 
Dirac point (i.e. n = everywhere in the sample) is 
doomed to fail no matter what. 

From our presented results in Sec. HI] (see also Tables U 
and ini , it is rather obvious that for the physics of in- 
trinsic plasmons to manifest itself experimentally, there 
would be a lower cut-off (T*) in the temperature be- 
low which (i.e. for T < T*) the electron- hole puddles 
would inhibit any observation of intrinsic plasmon behav- 
ior with the collective mode basically crossing over from 
being intrinsic for T >• T* to being extrinsic for T <C T* . 
Clearly, the temperature scale for this crossover is given 

by 

T* = Tpin*) ^ hvFVTm*/kB (46) 

where n* is the average puddle-induced carrier density 
in the system (i.e. n* is a function of rii). Thus, intrin- 
sic Dirac point behavior is essentially a high-temperature 
phenomenon with the intrinsic plasmon manifesting itself 
only for T ^ T* and then crossing over to extrinsic plas- 
mons of the puddle carriers for T <^ T*. The fact that 
the intrinsic Dirac point behavior is an effective "high- 
temperature" behavior has already been emphasized ear- 
lier in the context of approaching the Dirac point through 
transport measurements^i^, and in the current work, we 
establish the same qualitative finding for approaching the 
Dirac point through collective mode properties. 

In standard graphene on Si02 substrates, with sub- 
stantial impurity content in the environment, the typical 
puddle-induced carrier density n* ^ 10^^ cm~^, which 
gives T* ^ 1500 K, and obviously it makes no sense to 
discuss any experimental study of intrinsic plasmons in 
such "impure" graphene samples because the required 
temperature scale is impractically high. One can, of 
course, study extrinsic graphene in such samples by cre- 
ating a doped carrier density n > n* , and induced ex- 
trinsic graphene plasmons have been studied in such dis- 
ordered samples at high carrier density^'^'^. In sus- 
pended graphenei^^ or graphene on h-BN substratesi^, 
however, the environmental charged impurity density is 
very low (n^ < 10^ cm~^), and very low n* (^ 10^ — 10^ 
cm~^) has been reported with very sharp and narrow 
conductivity minimum at the charge neutrality point. 
Such ultrapure graphene samples (with typical mobil- 
ity ^Jim > 100,000 cm^/V-S) with very low puddle den- 
sity are the appropriate samples for studying intrinsic 
graphene. For n* = 10^(10^) cm'^, Tp ~ 50 K (15 K), 
and then the necessary condition (T ^ T*) for the man- 
ifestation of intrinsic plasmon collective behavior would 
necessitate T w 20 — 100 K, which is very reasonable 
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Material 



System 



High Temperature T - 



Comments 



Intrinsic 



7pi 



MLG 



Extrinsic 



Tp4 



'P' V 128(ln2)3Ti 



7pi 



(In 2)e-'q 
12refcHT 
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Intrinsic 



OP 
AP 



(qd)2 -ypj 



Double 
MLG 



Intrinsic- 
Extrinsic 



OP 
AP 



Two 

Extrinsic 



OP 

AP 



Same as results 
for double 
intrinsic layers 



For T -!> oo, 
chemical potential 
of MLG At ^ 



BLG 



Intrinsic 



\/2ujpi 



Extrinsic 



+ (41n2)T^ ^P^ 



V8 1 



8 In 2-2 Tp 
32 In 2 'P^ 
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Intrinsic 



OP 

AP 



2cVpi 



/2qdujpi 



VS{qd)2 -Ypi 



Double 
BLG 



Intrinsic- 
Extrinsic 



OP 
AP 



Two 

Extrinsic 



OP 
AP 



Same as results 
for double 
intrinsic layers 



For T -s> oo, 
Tp/T ^ 



2DEG 



Single layer 



/ 27vne~ q r 



3 qiikBT , 

4 Tve'^n J 



cxp 



-iksTq 



Double 
layer-'^ 



OP 
AP 



' 27re^ g(7il +712 ) 



/ ■7r(ni-t-n2) f ne^ \ 
V kgT \ nq ) 



■ n2 



"1+12 



/ TT / 27re dnin2 



exp[- 



l+exp[- 



ni+ri2 



TABLE I: Summary of high temperature analytical results. We denote ujpi = y ^ '-^^ — '■ — '-, "fpi 

only the leading order terms are kept in obtaining results for a double-layer 2DEG system. 



4(ln2)(fcaT)t;2q 



In 2 



Note 1: 



from an experimental perspective. We have already em- 
phasized in Section |lT] (and this is apparent in Tables U 
and ini that the extrinsic plasmon behaves precisely as 
an intrinsic plasmon (i.e. ujp cx y/T and a broadening 
7 oc l/y/T) for T ■> Tf, which then crosses over to 
the temperature-independent plasma frequency and ex- 
ponentially suppressed broadening, the characteristic fea- 
tures of extrinsic plasmon, at low temperatures T ^ Tp. 
This is what we would expect for high-quality graphene 
samples (n* ~ 10^ — 10^ cm~^) at the Dirac point with 
respect to their plasmon properties. The plasmon prop- 
erties will behave like those shown in Fig. [S] at higher 
temperatures (T » T* 10 — 50 K) to those shown in 
Fig.IDat lower temperatures (T < T* - 10-50 K). Such 
an observation will be a direct manifestation of intrinsic 
Dirac point behavior. 

Before concluding this section, we comment on the ex- 
pected plasmon level broadening (or damping) as man- 
ifested by the line width of the experimental plasmon 
peak. Our calculated damping (7) corresponds to the in- 
herent Landau damping induced level broadening which 
changes the plasmon peak from a pure delta-function like 
pole in the response function to an approximate broad- 



ened Lorentzian shape: 



6{U] — U!p) 



{uj - c^p)2 + r2 



(47) 



The broadening F in the Lorentzian plasmon peak in 
Eq. (|47p would, in general, have two contributions, aris- 
ing from Landau damping (7) and impurity broadening 
(7t) which can be written as 7* = S/(2t) where r is the 
transport relaxation time for impurity scattering as ex- 
tracted from the mobility. In the leading order, we can 
simply use: 



r = 7 -I- 7t 



(48) 



where 7 is the Landau damping calculated in SecHIlfsee 
e.g.. Tables HI and HI)) and 74 is the impurity-scattering 
induced level broadening which is given by: 



7* = fi/i^n) 



(49) 



For high mobility samples Tt oc is very large, and 
7f cx is small. Such a small impurity induced broad- 
ening will also be necessary for an unambiguous identifi- 
cation of the intrinsic plasmon — it is not enough to just 
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Material 


System 


Low Temperature T — > 


Comments 




7 




MLG 
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N/A 


N/A 
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V Y ' Tp-i 


'P»V (ln2)T 
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Two 
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OP 


/ 27re2g(ni+n2) 
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N/A 


7 exponentially 
suppressed 


AP 


V m/s:. Y ni+n2 


2DEG 


Single layer 




N/A 


7 exponentially 
suppressed 


Double 
layer 


OP 


/ 27re2g(ni+n2) 


N/A 


7 exponentially 
suppressed 


AP 


Y m/< Y '^l+"'2 



TABLE II: Summary of low temperature analytical results. We denote LOpxi — a/2 '' ■^^"^'^"^1 ^ ujpx2 ~ yi ^ ^Tm^tivpq ^ 



3 



/ 4(in2)(fcBT)e q ^ _ H. ^ Ji|. ^ £_2 ^ Note 1: results for a mixed intrinsic-extrinsic graphene double-layer system are 

obtained in the limit hvpq <^ ksT <^ Epi- For fcsT <^ hvpq, there is only an optical plamson mode ujop, determined by the 
extrinsic layer, in a mixed double-layer graphene system. 



have small values of only the intrinsic Landau damping 
7. Fortunately, the necessary condition for the suppres- 
sion of n* and T* (i.e. the suppression of puddles) also 
necessitates very large (small) values of Tt^-yt) since low 
impurity density (rii ~ n*) implies high mobility (since 
/i„i ~ l/?T-i) and low r*. Thus, very high-mobility sam- 
ples with very low puddle density, ensuring both 74 and 
n* to be small, would be necessary for the experimen- 
tal study of intrinsic plasmons. Fortunately, such high- 
quality samples already exist in the laboratory^, where 
the background impurity density (ui) is very low, thus 
ensuring that both T* and 74 are low enough for the ex- 
perimental investigation of intrinsic plasmons. 



IV. CONCLUSION 

We have provided a rather comprehensive theory for 
the intrinsic collective plasmon modes in graphene asso- 
ciated with the Dirac point within the dynamical finite- 
temperature random phase approximation. We have 
considered monolayer graphene, bilayer graphene and 
double-layer graphene (consisting of two MLG or BLG 
layers parallel to each other separated by a distance in 
the third direction). We have obtained within RPA ex- 
tensive analytical results for both plasmon dispersion and 
damping in the experimentally relevant long wavelength 
limit, and have provided detailed numerical results valid 
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for arbitrary wavevectors and temperatures. We have 
critically discussed the experimental feasibility for ob- 
serving the graphene intrinsic plasmon modes. 

Among our more important qualitative conclusions 
(arising from the theory developed in this paper) are: 
(i) although the intrinsic plasmon modes, being inher- 
ently finite-temperature excitations, are strongly Lan- 
dau damped even at low temperatures, they remain well- 
defined in the sense that the energy of the mode is larger 
than its damping in a very large regime of temperature, 
wavevector, and background dielectric constant; (ii) the 
intrinsic plasmon is inherently a "high-temperature" phe- 
nomenon since the Dirac point has no energy scale (and 
thus any finite temperature, no matter how low, is inher- 
ently a high temperature); (iii) closely connected with 
the last item is the corollary that the best experimental 
approach toward the observation of the intrinsic plas- 
mon is to study low-density extrinsic plasmon at tem- 
peratures high enough (i.e. T ^ Tp so that the effect 
of doping density is minimal); (iv) high-quality currently 
available graphene samples (either suspended graphene 
or graphene on h-BN substrates) with very high mobility 
should manifest clear-cut evidence for the intrinsic plas- 
mon (e.g. plasmon energy increasing as Vt and plasmon 
damping decreasing as VT with increasing temperature) 
if experiments are carried out at T sa 100 K with the gate 
voltage tuned to the nominal charge neutrality point; (v) 
the collective mode dispersion and damping for intrin- 
sic MLG and BLG plasmons is essentially identical (with 
Ldp cx \/qT and 7 cx q^/qT in both) except for numer- 
ical factors); (vi) our theoretically calculated analytical 
formula for plasmon dispersion and damping (both for 
intrinsic and for extrinsic plasmons) seem to agree with 
the full RPA numerical results essentially at all wavevec- 
tors and temperatures as long as our low-temperature 
analytical formula is used for T < Tp and the high- 
temperature analytical formula is used for T > Tp for 
extrinsic graphene; (vii) for double-layer systems, wc es- 
tablish that it should be experimentally possible to ob- 
serve both the optical and the acoustic intrinsic plasmon. 

Before concluding, we mention that the graphene plas- 
mon frequency is likely to be affected by interaction ef- 
fects even at long wavelength, unlike the long- wavelength 
plasma frequency in parabolic band systems which is 
protected by Galilean invariance and the associated /- 
sum rule so that only the band mass enters the def- 
inition of the long-wavelength plasma frequency, since 



graphene energy dispersion obeys Lorentz invariance. We 
believe that RPA is still an excellent approximation for 
the graphene plasmon properties (since RPA accounts for 
the long-range Coulomb potential correctly and nonper- 
turbatively) except that the velocity entering the expres- 
sion for the graphene plasmon mode should be modified 
to be the renormalized graphene velocity due to electron- 
electron interaction as calculated for example in Rcf. 
[9^. This simple modification in the spirit of Landau 
Fermi liquid theory should suffice to incorporate the lead- 
ing order interaction effect in the plasma frequency since 
the rs values characterizing Coulomb interaction strength 
in graphene are typically not too large. A more detailed 
theory for including interaction effects in the graphene 
plasmon properties is well beyond the scope of our work 
and would require substantial future theoretical efforts. 

We now conclude by discussing one important point 
which follows directly from our theoretical work with 
implications for graphene plasmonics. The recent in- 
terest in graphene plasmonics arises from the fact that 
graphene is a nanostructure allowing for very tight size 
confincmcntl^I and that the graphene (extrinsic) plas- 
mon has energy tunable by gate voltage (i.e. ujp ^ 
through the carrier density. In the case of 
intrinsic plasmon, however, the plasma energy goes as 
u}p ^ since basically the doping density n in the 

extrinsic formula gets replaced by the thermal electron- 
hole excitation density (i.e. n ^ T^). This implies that it 
should be relatively easy to tune the intrinsic plasmon en- 
ergy simply by changing temperature while sitting at the 
charge neutrality point,. Since the temperature depen- 
dence (~ VT) of the intrinsic plasmon frequency is much 
stronger than the doping density dependence (i.e. n^/"*) 
of the extrinsic plasmon frequency, it may be more conve- 
nient to use the thermal tuning of the plasma frequency 
than the gate voltage tuning considered so far. In addi- 
tion, the intrinsic plasmon should be a relatively strong 
and well-defined mode at room temperatures (T ^ 300 
K) in high-mobility graphene samples, making it an in- 
teresting candidate for possible plasmonic applications. 
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